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Abstract 

We consider evaporation of pure liquid drops on a thermally conductive substrate. Two evapo- 
rative models are considered: one that concentrates on the liquid phase in determining evaporative 
flux, and the other one that centers on the gas/vapor phase. A single governing equation for the 
evolution of drop thickness, including both models, is developed. Experiments are used to estimate 
relevant parameters. We show how the derived governing equation can be used to predict which 
evaporation model is appropriate under considered experimental conditions. 

PACS numbers: 68.03.Fg,68.15.+e,47.55.Ca,47.55.np 
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Evaporating thin films and drops are present in numerous natural situations and ap- 
plications of technical importance. Coated liquid films, for example, are often left to dry 
by evaporation. The residual films, whose thickness may vary from millimetric in the case 
of paints to nanometric for photoresist films in semiconductor applications, are often de- 
sired in a uniform state. However, various kinds of instabilities, many of them driven by 
evaporation-related mechanisms, often occur. Evaporative sessile drops are perhaps even 
more interesting since nonuniform drop thickness and presence of contact lines (separating 
liquid, gas, and solid phase) lead to new effects. One of these is a possibility of nonuniform 
evaporation along the gas-liquid interface, leading to temperature gradients and related 
Marangoni effects. These effects are crucial in a number of problems, including the so-called 
coffee-stain phenomenon involving deposition of solid particles dissolved in the liquid close 
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to a contact line [I], and its' numerous applications, such as analysis of DNA microarrays [2j. 

Despite its apparent simplicity, the problem of an evaporating drop on a thermally con- 
ducting solid substrate involves a number of physical processes, including mass and energy 
transfer between the three phases, diffusion and/or convection of vapor in the gas phase, 
coupled to the complex physics in the vicinity of a contact line. So-called '2-sided' models 
include the processes both in the liquid and in the gas phase, but lead to a mathematical 
formulation of significant complexity, even when the solid phase and contact line issues are 
not considered |3| . Therefore, most of the researchers in the field have chosen to use various 
simplifications, allowing them to reduce the complete problem to a tractable mathematical 
formulation. As we discuss below, these simplifications are based on estimating the impor- 
tance of various physical processes, and lead eventually to models which concentrate on the 
processes in one of the phases (gas or liquid). The estimates, however, involve quantities 
which are either not known, or not known precisely enough. We will show that various as- 
sumptions may lead to models which can produce qualitatively different and contradictory 
results. The answer to the question 'How do drops evaporate?' is still not known. 

The complete 2-sided model can be simplified by realizing that thermal conductivity and 
viscosity of vapor are small compared to the ones of the liquid. In addition, assuming that 
the gas phase is convection-free, one reduces the 2-sided model to the so-called '1.5-sided' 
model which includes the processes in the liquid, and the diffusion of vapor in the surrounding 
air 4j. An estimate of a typical diffusion time scale = l 2 /D, involves the relevant thickness 
of the gas phase, I, and the diffusion constant for e.g. water vapor D ~ 10~ 5 m 2 / s. Assuming 
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for a moment that I is on millimeter scale (comparable to a typical thickness of a drop), leads 
to td pd 10 -2 s. The argument that is much shorter than a typical time scale involved in 
drop evolution has been used to reduce the diffusion equation for vapor concentration, c, to 



a Laplace equation [5]. Furthermore, assuming that evaporation process itself is extremely 
fast [6] allows to completely ignore the processes in the liquid for the purpose of finding 
the mass flux. Therefore, the problem is simplified to V 2 c = in the gas/vapor domain. 
Concentrating now on the part of domain close to the contact line, one realizes electrostatic 
analogy of finding an electric field (mass flux, J) in the vicinity of a 'lens' shaped conductor 
(the drop), where c plays the role of electrostatic potential [1]. Typically this lens model 
assumes a pinned (stationary) contact line, although this assumption is not crucial since 
typically the contact line, even if mobile, evolves slowly Q. One important outcome of the 
model is that, under some additional simplifications, J oc l/h x , where h is the drop thickness, 
and A = A(9), being the contact angle [8|. Clearly, J diverges at the contact line, and 
various (often ad-hoc) procedures have been used to regularize the problem. An extensive 
modeling using lens- type model has been carried out, implementing both finite-element and 



lubrication- type approaches 



9fl . Various versions of this model have been used to predict 



evaporative behavior of alcane [lOj] and colloidal drops [8|, LL1[ and the temperature along 
liquid-gas interface among others. 

Turning now to the liquid phase, one realizes that the evaporation is limited by two 
physical processes: heat diffusion through the liquid supplying heat to the interface, and 
evaporation itself. In a simple model |3(, these two processes can be related via Biot's 



number, Bi = Kp T Ld /(p v k), where K = ap v {T { ) / ^2TiRgT sat and p T = Lp sat /(R g T^ at ). 
Here, p sa t is the saturation pressure, R g is the universal gas constant divided by the molar 
mass, L is the latent heat of vaporization, do is drop thickness, p v is vapor density, k is 
liquid heat conductivity, and T sat are the interface and saturation temperatures, and a 
is the accommodation coefficient, describing probability of phase change. The limit Bi — > 
implies that the temperature of the liquid/gas interface tends to the temperature of the 
solid, evaporation itself proceeds in reaction-limited regime, and the interface is in a state 
of non-equilibrium Q] . — > oo, on the other hand, indicates that the evaporation is much 
faster than the heat diffusion through the liquid, and evaporation proceeds in liquid heat 
diffusion-limited regime, with the interface being in an equilibrium state. 

While most of the quantities entering the definition of Bi are well known, the value of 



a is questionable. A variety of a's in the range 0(1) — O(10 -6 ) have been used in the 
literature, often without much justification. We note that the theoretical predictions that 
a G [10~ 2 , 1] for water have been found to grossly overestimate the volatility, with values 
in the range [10 -6 , 10 -4 ] being more realistic 13j]. In our experiments described below, 
we find a ~ 10~ 4 , consistently with other experimental results involving drops exposed to 
ope, atm0S p h ,e Q. T, li8 a gives Bi w 10 2 , suggesting that evaporation proceeds in 
the reaction-limited regime. Further insight can be reached by considering Bi as the ratio 
of the relevant time scales involved in heat diffusion in the liquid, ti, and the one related 
to the evaporative process itself, t e . Using ti = do 2 / 'k, where k is thermal diffusivity of 
liquid, and do = 0.5 mm, one finds t\ ~ lsec, giving t e ~ 10 2 sec. Similar value is obtained 



by using t e = (2pdoL)/(kAT), where AT is the appropriate temperature scale |14J ] . Since 
t e ^> ti ^> td, one may consider a model where relevant limiting mechanism is the evaporative 
process itself, and not the diffusion of vapor in the gas. This 'non-equilibrium one sided' 
(NEOS) model has been extensively used for a variety of problems involving evaporative thin 



films 



151 ]. but only few works have applied it to evaporating drop problem [111 . 
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note that the argument outlined here for use of the NEOS model is based on the assumption 
of relatively small relevant thickness, I, governing diffusion in the gas phase. Larger I would 
lead to larger td, and thus it would be unclear which model is more appropriate. We note 
that, if td is large, the assumption of steady-state formulation for the gas concentration is 
questionable. 

In this letter we show computationally that these two commonly used models for evap- 
oration may lead to inconsistent results regarding volume loss, time evolution of the drop 
size, and in particular regarding the temperature gradient along the liquid-gas interface. 
We show that each model requires the use of unknown quantities, making it difficult to 
decide on which one is more appropriate. Therefore, we carry out experiments to estimate 
these quantities and compare the two models, with the goal of reaching a conclusion regard- 
ing an appropriate one, at least for the experiments considered. To our knowledge, this is 
the first time these two models have been compared directly against each other using data 
appropriate to a particular experiment. 

The mathematical model consists essentially of Navier-Stokes equations for the liquid 
phase coupled with the energy equation for the solid, and appropriate expression for J, in 
the spirit of the models reviewed in 15[]. We use lubrication approximation although for 
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some considered problems (such as a dionized water (DIW) drop on a silicon (Si) substrate) 
the contact angle is relatively large (0 ~ 40°). This approach is supported by finite- 
element simulations that show that even for large contact angles considered here, lubrication 
approach leads to reliable results [5]. Solid-liquid interaction is included using disjoining 
pressure model with both attractive and repulsive van der Waals (vdW) terms, leading to 
a stable equilibrium thickness, which can be related to commonly used precursor film. In 
cylindrical coordinates and assuming azimuthal symmetry, we obtain the following 4th order 
PDE for the film thickness, h(r, t) 
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which is put in nondimensional form by choosing do and d\jv as the length and time scale, 
respectively [v is the kinematic viscosity). We used a similar formulation to consider insta- 
bilities of evaporative isopropyl alcohol (IPA) drops [ijj]. The main difference of the present 
case (in addition to the geometry) is keeping the evaporative flux J(h) explicitly in the 
formulation, so that Eq. ([I]) can be used for any evaporative model. Eq. ([T]) includes the ef- 
fects due to viscosity, evaporation, capillarity, vapor recoil, Marangoni and vdW forces, and 
gravity. All nondimensional quantities are defined in [14j |. and here we just enumerate the 
most important ones: E, S, M, P, A, G correspond to evaporation number, nondimensional 
surface tension, Marangoni number, Prandtl number, nondimensional Hamaker constant, 



and gravity. The material parameters are listed in 



Next we concentrate on the key point, and that is the evaporative flux, J. For the lens 
model, we use J(h) = x/h x , where the exponent A can be approximated by A = 0.5 — Q/tt j^]. 
For the NEOS model, we use J(h) = l/(h + W + /C) [14], where W accounts for finite 
thermal conductivity of the solid substrate and K, = Bi~ l . For the case of DIW, /C is 
typically large, as discussed before, and therefore J(h) only weakly depends on h. The 
crucial parameters, volatility \ an d a are obtained from the experiments described next. 

The experiments are carried out at room temperature and in open atmosphere, using a 
goniometer (KSV CAM 200) which consists of a camera, light source, and static deposition 



platform. Pure liquids and silicon wafers of semiconductor grade smoothness are used. 
Figure [T^) shows an example of an image of evaporating DIW drop. The build-in software 
analyzes the images and yields the height and the radius of a drop. Figure [Tb) shows the 
resulting radius and volume of evaporating drop as a function of time. In agreement with 
other works, we find linear decrease of the volume for the considered time interval 
We note that we do not observe either significant stick-slip motion, nor contact line pinning 
that is often assumed to occur jj], 9]. Possibly, these effects are not seen due to small surface 
roughness of the Si wafers and due to purity of the liquids. 
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FIG. 1: (a) Snapshot of evaporating DIW drop on Si substrate; (b) Volume and radius of the drop 
as a function of time (experimental results). 



In order to extract \ an d ot from the experimental data, we proceed as follows. The 
considered experimental time frame (typically the first 106 sec), is split into subintervals, 
and then the mass loss in each subinterval is used to obtain appropriate values of x an d ot 
by integra ting J over the drop surface (more details regarding this procedure will be given 
elsewhere l_7j] ) . We expect that this method is more accurate than simply using the dry 
out time to estimate the flux, since at late stages of evolution, the evaporation rate may be 



reduced 



181 ] . Figure [2] shows the results. For the purpose of verification, we have applied the 



51 and obtained similar values for 



same approach to the already available experimental data 
a and \. With these values available, we are now almost ready to proceed with numerical 
simulations of the drop evolution governed by Eq. The only missing ingredient is the 
precursor film thickness, b. Although it is typically sufficient to use 6 C 1, here we are 
more careful since quantitative agreement is desired. We are governed by the requirement 
that r(t) and V(t) should not depend on b in any significant manner, and have found 
that for d Q b < 0.625 fim this requirement is satisfied. Coincidentally or not, this value is 
consistent with the equilibrium adsorbed film thickness dob e for which evaporation stops due 
to attracting vdW forces [if]]. For the appropriate parameters, we find d§b e ~ 0.5 j[/m. We 
note that consideration of vdW forces automatically regularizes otherwise singular expression 



for J in the lens model [19|. Without inclusion of vdW effects, an additional externally added 
regularization of J has to be included in order to correctly compute the total mass flux. 
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FIG. 2: Volatility for (a) lens and (b) NEOS models for DlW/Si, leading to % ~ 5.5 • 1(T 2 and 
a ps 3.4 • I(T 4 (lines). 

The numerical simulation of Eq. (pO) is carried out using second order accurate implicit 
scheme which is an extension of the one used in [fjj] to cylindrical geometry. All simulations 
use as an initial condition the steady-state solution of Eq. ([1]) obtained by removing the 
evaporative terms; in experiments, this configuration is reached after a very short time. 
Figure Eh) shows the volume of the evaporating drop as a function of time. We find very 
good agreement between the experimental results and the NEOS model, while the lens 
model overestimates the volume loss due to evaporation. Therefore we conclude that, at 
least for the problem considered, the NEOS model predicts better the volume loss compared 
to the lens model. We discuss this difference, as well as additional experiments/simulations 
involving evaporation on different solid surfaces and under modified experimental conditions 
elsewhere [17]. 
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FIG. 3: (a) Volume of an evaporating DIW drop: comparison of the models with the experiment, 
(b) Temperature of the liquid-gas interface predicted by the two models at the final time shown in 
(a). 

Figure [3b) shows that the two models predict qualitatively different temperature profiles 
along the liquid-gas interface. An increase of temperature as one moves from the center 



in the NEOS model is the consequence of the fact that the heat supplied from the solid is 
larger than the heat loss due to evaporation. The lens model, on the other hand, predicts 
significantly larger evaporative flux in the contact line region, therefore leading to a sharp 
decrease of temperature. An increase of temperature as one moves away from the drop center 
is consistent with the results [3] obtained using the lens model, and similar values of 
and thermal conductivities of the liquid and solid, although under the assumption of a 
pinned contact line. However, the results presented here predict a 'stagnation point', where 
the temperature gradient changes sign. Presence of a stagnation point, although based on 
different physical grounds, was discussed earlier as one of the necessary ingredients required 
for formation of ring-like deposits occurring during evaporation of colloidal drops [8]. We 
note that simulations assuming the lens model for smaller values of G lead to monotonically 
decreasing temperatures along the liquid-gas interface as one moves away from the drop 
center, in full agreement with the earlier results fl2. I^ol]. 

Next, we apply these two models to a more volatile IPA drop. This configuration is 
characterized by ~ 6° leading to very thin drops. Due to this fact, in our experiments we 
cannot follow accurately the evolution in the manner it was done for DIW drops. Therefore, 
we use the dryout time, which we can estimate with a reasonable accuracy, to obtain the 
required values of a and x- F° r V(Q) ~ 3.2 mm 3 , we find td ry ~ 135 sec. Applying the two 
evaporative models leads to a = 9.3 • 10~ 4 and x — 5.7 • 10~ 3 . We note that similar values of 
X were obtained in experiments with pure alcanes [lOj • Implementing; these values, we solve 
Eq. ([1]) using the parameters appropriate for IPA and Si [14j . 

Figure 0k.) show the resulting evolution of the drop radius, R(t), and the temperature 
along the liquid-gas interface (b). First, we notice dramatically different evolution of R(t) 
for the two models. Considering the temperature profiles in 0b) provides immediate under- 
standing of this difference. The Marangoni forces act in the opposing directions, leading 
to a very different evolution. For example, for the lens model, the Marangoni forces act 
outwards, leading to initial increase of the drop radius (shown in Fig. 0k) despite the loss of 
mass due to evaporation. The difference between the temperature profiles predicted by the 
two models is much more significant in the case of IPA compared to DIW due to its larger 
volatility. We discuss the details of model predictions regarding evolution of IPA drops in 
more details elsewhere ~v\ . 

We find it intriguing that the two commonly used models for evaporation produce results 
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that are qualitatively different, regarding the evolution of drops' volume, radius, and, in 
particular, the liquid-gas interfacial temperature. The resulting Marangoni forces may act 
in the opposing directions. Therefore, it appears that the arguments involving, e.g., the 
influence of Marangoni forces on formation of particle deposits next to the contact line have 
to be carefully reconsidered. We hope that the results presented here will encourage more 
elaborate experiments, possibly involving direct measurement of the interfacial liquid-gas 
temperature, an ultimate test for any evaporation model. 
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FIG. 4: IPA drop: (a) Radius R(t) as a result of the two models; (b) Temperature profiles of the 
liquid-gas interface at the final time shown in (a). 
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